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Abstract 

The dynamic rheology of a polymer melt composed of short chains with ten beads between 
rapidly oscillating plates is investigated for various oscillation frequencies by using the hybrid 
simulation of the molecular dynamics and computational fluid dynamics. In the quiescent state, 
the melt is in a supercooled state, and the stress relaxation function G{t) exhibits a stretched 
exponential relaxation on the time scale of the a relaxation time Tq (the structural relaxation 
of beads) and then follows the Rouse relaxation function characterized by the Rouse relaxation 
time tr (the conformational relaxation of polymer chains). In the rapidly oscillating plates, non- 
uniform boundary layer flows are generated over the plate due to inertia of the fluid, and the local 
rheological properties of the melt are spatially varied according to the local flow fields. The local 
strain and local strain rate of the melt monotonically decrease with the distance from the plate at 
each oscillation frequency of the plate, but their dependencies on the oscillation frequency at a fixed 
distance from the plate vary with the distance. Far from the plate, the local strain decreases as 
the oscillation frequency increases such that the dynamic rheology deviates from the linear moduli 
at the low oscillation frequencies rather than high oscillation frequencies. On the contrary, near 
the plate, the local strain rate increases with the oscillation frequency such that the shear thinning 
is enhanced at high oscillation frequencies. In close vicinity to the plate, the dynamic viscosity is 
mostly independent of the oscillation frequency, and the shear thinning behavior becomes similar 
to that observed in steady shear flows. We show the diagram of the loss tangent of the melt for 
different oscillation frequencies and local strain rates. It is seen that the melt generates three 
different rheological regimes, i.e., the viscous fluid regime, liquid- like viscoelastic regime, and solid- 
like viscoelastic regime, according to the oscillation frequency and local strain rate. Non-linear 
rheological properties are also investigated by the spectrum analysis and the Lissajous-Bowditch 
curve. It is found that the fractional amplitude of the higher harmonics to the linear harmonics is 
suppressed within the boundary layer due to the non-slip boundary on the oscillating plate. We 
also find that the melt exhibits inter-cycle shear thinning between different positions but exhibits 
intra-cycle shear thickening at a fixed position in the vicinity of the plate. 

PACS numbers: 

Keywords: dynamic rheology, polymeric fluid, stokes layer, multiscale modeling 
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I. INTRODUCTION 



Glassy polymeric fluids have complicated shear-dependent dynamic rheology. In steady 
shear flows, glassy fluids are highly viscous, and shear thinning occurs at sufficiently large 
shear rates due to chain elongations.-i^ In unsteady flows, glassy fluids exhibit elastic be- 
haviors if the characteristic time scale of polymer dynamics is comparable to or larger than 
that in the flow system. The viscoelastic property can be measured by the shear moduli, 
i.e., the storage modulus G' for the elasticity and the loss modulus G" for the viscosity. 
The shear moduli can be measured under uniform oscillatory shear flow with a finite shear 
strain, both experimentally or numerically. In general non-uniform flows, however, the local 
rheological properties become heterogeneous depending on the local flow variables. Thus, 
the rheological behaviors of glassy polymeric fluids in highly non-uniform flows are so com- 
plicated that the theoretical or experimental approaches to this problem are very difficult. 
It is also difficult to predict the flow behaviors of such fluids because the reliable constitutive 
equations are not known in general, although there is an important accumulation of both 
experimental and theoretical works to construct them.-^^ 

In the present paper, we investigate the dynamic rheology of a model polymer melt 
composed of short chains between rapidly oscillating plates by using the hybrid simulation 
of molecular dynamics (MD) and computational fluid dynamics (CFD). The temperature of 
the melt is so low that the glassy behavior is observed in the stress relaxation function in the 
quiescent state. Non-uniform oscillatory shear flows are generated over the oscillating plate 
at sufficiently large oscillation frequencies due to inertia of the fluid via the term pdv/dt, 
which is sometimes called the transient force. Thus, heterogeneous rheological behaviors 
arise according to local flow variables. 

Oscillatory shear flows under the transient force for the viscoelastic fluid have been in- 
vestigated by several researchers so far. Scharg analytically solved the flows of a linear 
viscoelastic fluid^. Dunwoody investigated the inertia effect for a weak non-linear viscoelas- 
tic fluid with the K-BKZ model by using a perturbation analysis^. Yosick et al. performed 
numerical analysis for a non-linear viscoelastic fluid with the Berkeley model. Ding et al. 
investigated the viscous dissipation under the transient force in the temperature field for 
linear viscoelastic fluids*. In previous works, some constitutive model equations have been 
used to calculate flow profiles. In the present hybrid simulation, no constitutive model is 
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required to obtain a local stress. The local stress is generated by a local MD simulation 
according to the local flow variables. 

In the hybrid simulation method, the macroscopic flows of the melt are calculated by using 
the CFD scheme; however, instead of using any constitutive equations, the local stresses of 
the melt are calculated by using molecular dynamic simulations of polymer chains according 
to the local strain rates.-""— The basic idea of the present hybrid simulation method was 
first proposed by E and Engquislr^'i^, where the heterogeneous multiscale method (HMM) is 
presented as a general methodology for the efficient numerical computation of problems with 
multiscale characteristics. The HMM has also been applied to the simulation of complex 
fluids.— Equation-free multiscale computation was also proposed by Kevrekidis et al. on 
the basis of a similar idea and has been applied to various problems.— De et al. have 
developed a hybrid method, called the scale bridging method, which can correctly reproduce 
the memory effect of a polymeric liquid and demonstrated non-linear viscoelastic behavior of 
a polymeric liquid between oscillating plates.— The multiscale simulation based on a similar 
idea has been also applied to rarefied gas flows, recently.— 

In what follows, we characterize the linear dynamic rheology of a model polymer melt by 
using the MD simulation in Sec. II. Then, we investigate the dynamic rheology of the melt 
in non-uniform oscillatory shear flows under the transient force between rapidly oscillating 
plates in Sec. Ill, where the simulation method, velocity profiles, shear moduli, and LAOS 
analysis are introduced. Finally, we summarize the results in Sec. IV. 



II. LINEAR DYNAMIC RHEOLOGY OF A MODEL POLYMER MELT 

We consider a model polymer melt composed of short chains with ten beads of a uniform 
density po ^^nd a uniform temperature Tq. The number of bead particles on each chain is 
represented by A^b- Thus, N\j = 10. All of the bead particles interact with a truncated 
Lennard- Jones potential defined by^, 

j Ae[{a/ry'-{a/r)^]+e (r < 2^^), 
I (r > 2i/V). 

By using only the repulsive part of the Lennard- Jones potential, we may prevent spatial 
overlap of the particles. Consecutive beads on each chain are connected by an anharmonic 
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spring potential, 

U^{r) = -h,Rl\n[l~{r/R,f], (2) 

where A;c=30e/(7^ and i?o=l-5cr. The temperature of the melt is A;BTo/e=0.2, where ks is the 
Boltzmann constant. The number density of the bead particles is pQ/m=l/(y'^, where m is 
the mass of the bead particle. With this number density and temperature, the configuration 
of the bead particles becomes severely jammed, resulting in a complicated non- Newtonian 
viscosity and long-time relaxation phenomena characteristic of glassy polymers.->^ Here- 
after, unless otherwise stated, we measure the physical quantities with units of length a, 
time \fma^e, and temperature e/fc^. 

In this section, we clarify the linear viscoelastic properties of the model polymer melt. 
Figure [1] shows the stress relaxation function Gii) in the quiescent state. The stress re- 
laxation function Gif) of the model polymer melt was obtained in Ref. |20| (although we 
recomputed in the present study in order to calculate the accurate linear complex 
moduli as seen below). The details of the molecular dynamic simulation can be found in 



Ref. 



2d . The stress relaxation function Gif) is calculated as 



G{t) =< U,y{t + to)n.,(to) > /ksTV, (3) 

where U^y is the space integral of the microscopic stress tensor in the volume V. In 
the figure, the a relaxation time Tq, and the Rouse relaxation time r^j in the quiescent 
state, where r„=310 and r/j=6 x 10^, are also plotted. The a relaxation time Tq, represents 
the characteristic time of the structural relaxation of bead particles and is calculated from 
the van Hove self-correlation function, and the Rouse relaxation time tr represents the 
characteristic time of the conformational relaxation of polymer chains and is calculated 
from the time- correlation function of the end-to-end vector of each polymer chain. After 
the early oscillating behavior for t < 10, which corresponds to the vibrations of the bonds 
of bead particles on each chain, G{t) exhibits the stretched exponential relaxation of the 
Kohhausch- Williams- Watts (KWW) form G,(t), 

a(t)=G'oexp[-(t/r,)^]. (4) 

G{t) can be nicely fitted to Eq. (jlj) with Go=5, c=0.5, and Ts = 90 0.33rQ,) for 
1~ ^ ~ 10^- Ultimately, G{t) follows the Rouse dynamics characterized by the Rouse 
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relaxation time tr. The Rouse relaxation function Gii{t) is written as 



PqT 

iVb 



exp(-2t/rp), 



(5) 



p=i 



where Tp/rji = sin^(7r/2A^b)/ sin^(vrp/2A^b) for p = 1, ■ ■ ■ ,A'b. Thus, in the quiescent 
state, the stress relaxation function G{t) can be approximately described by the sum of 
the stretched relaxation function Gg in Eq. and the Rouse relaxation function G^lt) 
in Eq. ((51), i.e., G(t) ~ Gs(t) + Guit), except the period of initial oscillating modes, and 
decreases so rapidly that it becomes negligible in the late stage, while G_R(t) is so small as 
to be negligible in the early stage but can describe the late stage of G{t). 

The frequency- dependent shear moduli, i.e., the storage modulus G'^uj) and the loss 
modulus G"{u)), in the linear response regime are obtained by the Fourier transform of the 
stress relaxation function G{t) in the quiescent state as 



G'{uj) and G"{u) represent the elasticity and viscosity of the melt, respectively. Figure 
[2] shows the linear storage modulus G"(w) and loss modulus G"{ijj) versus the frequency u. 
The crossover of G'{u) and G"{uj) is observed aX u t^^. The model polymer melt is rather 
elastic {G' > G") for u > r~^, while it is rather viscous {G" > G') for u < r~^. At the low 
frequency, i.e., u < r^^, the loss modulus G" is quite dominant and almost proportional to to. 
This indicates that, as the frequency to is smaller than r^^, the dynamic viscosity rj' defined 
as 7]' = G"/u) is approximately constant so as to be a simple viscous fluid with a constant 
viscosity rj'. The shear moduli calculated from the sum of the stretched exponential form 
Gs{t) and the Rouse relaxation function Gpi{t) are also plotted in the figure. It is seen that, 
in the linear response regime, the moduli of the model polymer melt can be well described 
by those calculated from the sum of Ggif) and GrIJ:), except for the high frequency regime 
as w ^ T^^, which corresponds to the oscillation mode of G(t) in the early stage (See Fig. 




oo 



G{t) sinut dt, 



(6a) 




(6b) 
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III. THE EFFECT OF INERTIA OF THE FLUID 



A. Problem and simulation method 

In the rapidly oscillating plates, the oscillatory shear flow becomes non-uniform due 
to inertia of the fluid via the term pdvx/dt, say the transient force, such that the local 
rheological property spatially varies according to the local velocity fleld. In this section, 
we investigate the dynamic rheology of the melt in the non-uniform oscillatory shear flow 
between rapidly oscillating plates (see Fig. [3]^a)). The upper- and lower- plate start to 
oscillate in the x-direction at a time t=0 as, respectively, 

Vw(t) = T'^oHuocos^uot), (7) 

where Fq is the amplitude of the strain of the system and H is the half of width between 
the upper- and lower- plate. We assume that the macroscopic quantities are uniform in 
the X- and ^-directions, d/dx=d/dz=0. The macroscopic velocity Va is described by the 
following equations, 

dvx da^y 

and Vy=Vz=0, where t is the time and a^y is the shear stress. We also assume a non-slip 
boundary condition on each plate. In the present study, we flx the width between the plates 
at 2if=5000. For the strain amplitude of the system Gq, Fo=0.5 is the value that is mainly 
utilized, while Go=0.02 is a subsidiary used for the comparison. However, the oscillation 
frequency of the plate Uq is widely varied in order to investigate the effect of the changing 
oscillation frequency on the rheological property of the melt. The oscillation frequencies uo 
used in the present simulations are shown by diamonds O on the upper-axis around uq ~ 
in Fig. H 

We solve this problem by using a hybrid simulation of the molecular dynamics (MD) 
and computational fluid dynamics (CFD). The details of the method can be found in Ref. 
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Here, we briefly explain the hybrid simulation method. We calculate the macroscopic 
velocity in Eq. (|8]) by using a usual flnite volume scheme with a uniform mesh system 
(see Fig. Mjo) )• However, instead of using any constitutive equation for a^y, we calculate 
the local stresses in small MD cells associated with each mesh interval according to the 
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local strain rates, which are calculated at the CFD level, at each time step of the CFD 
simulation.— 1^ The MD simulations are performed in each MD cell for the duration of the 
time-step size of the CFD calculation, and the molecular configurations obtained in each 
MD cell after one MD run are memorized as the initial configurations of molecules for the 
next MD run in each MD cell (see Fig. ^c) ). By using this method, one can reproduce 
the memory effect caused by the slow dynamics of the molecular conformation correctly. In 
the present study, we divide the lower half between the plates into 128 mesh intervals with 
a mesh size of Ax=19.5 and use the symmetric condition at the middle between the plates 
for the CFD calculation. For the MD simulation, we use a cubic MD cell with a side length 
^MD=10 so that each MD cell contains 1000 bead particles. The ratio of the mesh size of 
CFD Aa; to the side length of the MD cell Imb, which represents the efficiency of the hybrid 
simulation compared to the full MD simulation, is Aa;/ZMD=l-95. The time-step size of the 
CFD calculation At and that of the MD calculationAr are fixed at At=l and Ar=0.001, 
respectively. Thus, one thousand time steps are carried out in each MD run at each time 
step of the CFD simulation (M = 1000 in Fig. [3] (c)). 

B. Velocity profile 

FigureH] (a) shows snapshots of the velocity profile of the melt for a;o=0.0015 and ro=0.5. 
Due to the transient force, the amplitude of the oscillatory velocity rapidly decreases with 
distance from the oscillating plate, and a thin boundary layer forms over the oscillating 
plate. The thickness of the boundary layer 1^, which is defined as v^iy = h)/vo = e~^, is 
also shown at different oscillation frequencies uq for ro=0.02 and 0.5 in Fig. H] (b). The 
solid line shows the thickness of the boundary layer for the linear viscoelastic fiuid with a 
linear storage and loss modulus shown in Fig. O As is seen in Fig. |2l at the low frequency 
(i.e., CO < T^^), the storage modulus G' is much smaller than the loss modulus G" and the 
dynamic viscosity 1]', which is calculated as G" /u, becomes almost constant. Thus, the slope 
of w"^/^ for the Newtonian fiuid arises at the low frequency in Fig. H] (b). The thickness 
of the boundary layer of the present polymer melt is close to that of the linear viscoelastic 
fiuid at low oscillation frequencies for the small strain amplitude of the system, Fq = 0.02; 
however, as the oscillation frequency increases, the boundary layer of the melt becomes much 
thinner than that of the viscoelastic fiuid. The thickness of boundary layer for the large 
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strain amplitude, Fq = 0.5, is thinner than that for the small strain amplitude, Tq = 0.02 
at any of the oscillation frequencies. These features are caused by the shear thinning of 
the melt because the local strain becomes larger near the oscillating plate as the oscillation 
frequency Uq and the strain amplitude of the system Tq increase. Due to the emergence of 
the boundary layer, the rheological properties of the melt become more complex in nature, 
as we see below. 

C. Storage and loss modulus 

In this subsection, we investigate the "local" rheological properties of the melt in the slab 
in terms of the storage modulus G' and loss modulus G". The local moduli are calculated 
from the first harmonics of the Fourier series of the time evolutions of the local shear stress 
and strain. By using the Fourier transform of the time evolution of the local strain 7(2/, t) 
and selecting the mode of the oscillation frequency of the plate coq, we can approximate the 
time evolution of the local strain 7 in the form of 7 = 70 (y) cos((X'ot + '?/'(?/)). Here, 70 (y) is 
the strain amplitude, and ip{y) is the phase retardation. In the same way, the local shear 
stress can also be written clS (J xy — '{y) cos{uJot + — o" sin(a;o^ + i')- The local storage 
modulus and loss modulus are obtained as G'{y) = o-' (y) / ■yo{y) and G"{y) = o-" (y) / Joiy) , 
respectively. We note that the contribution of the higher harmonics is also important in 
the large amplitude oscillatory shear (LAOS) regime and is actually detected in the present 
simulations. However, the amplitude of the higher harmonics is smaller than that of the first 
harmonics; the fractional power spectrum of the higher harmonics of the local shear stress 
is at most 2.8% in the present simulations. Thus, the storage and loss modulus represent 
the basic viscoelastic properties of the melt in the slab. (The non-linear effect of higher 
harmonics is investigated in the next subsection.) In the present subsection, we fix the 
strain amplitude of the system as Fq = 0.5 and change the oscillation frequency of the plate 
Uo variously in order to investigate the dependency of the local rheological property on the 
oscillation frequency of the plate. 

Figure [5] shows the dependency of the storage and loss modulus on the oscillation fre- 
quency, G' and G" vs. ujq and that of the strain and strain rate amplitudes on the oscillation 
frequency, 70 and 70 vs. Uq, at a fixed position. Figure (a) shows the storage and loss 
modulus versus oscillation frequency and the strain and strain rate amplitudes versus os- 
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cillation frequency at far distances from the plate and (b) those at the near distances from 
the plate. It can be seen that, at the high oscillation frequency, say uo = 0.01, the stor- 
age modulus G' is larger than the loss modulus G" at large distances from the plate, e.g., 
y = 1309 (Blue) and 771 (Green) in Fig. [S] (a), while G' is smaller than G" as the distance 
decreases, e.g., y = 381 (Red) in Fig. \5\ (a) and at smaller distances in Fig. \5\ (b). In close 
vicinity to the plate, e.g., ?/ = 10 in Fig. |5] (b), the viscosity is quite dominant, G" ^ G'. 
Thus, the local rheological properties of the melt vary considerably according to the local 
flow fields such that three different rheological regimes, i.e., the viscous fluid regime, the 
liquid-like viscoelastic regime, and the solid-like viscoelastic regime, are formed over the 
rapidly oscillating plate. On the contrary, at a low oscillation frequency, say = 0.001, the 
differences of the storage modulus G' and loss modulus G" are not as large as those at a high 
oscillation frequency at any distance, and the loss modulus G" is larger than the storage 
modulus G'. The differences of local modulus between the different distances are also not as 
large as those at a high oscillation frequency. Thus, the local rheological properties of the 
melt vary moderately between the plates, and the melt has liquid-like viscoelastic behaviors 
throughout the slab. 

At large distances from the plate [Fig. [5] (a)], the local strain 70 monotonically decreases 
as the oscillation frequency ojq increases because, as we have seen in Fig. HJ a thin boundary 
layer forms over the rapidly oscillating plate, and the thickness of the boundary layer becomes 
thinner as the oscillation frequency uq increases. Thus, the local moduli deviate from the 
linear moduli more at the lower oscillation frequencies due to the shear thinning, while they 
are rather close to the linear values at the high oscillation frequencies. 

The shear thinning behaviors of the local moduli to the local strain amplitude 70 are 
shown in Fig. [6l It is seen that the storage modulus G' decreases when the local strain 
7o is larger than about 2%, while the loss modulus G" starts to decrease at a larger strain 
amplitude (70 > 2%). The non-monotonic behavior of the loss modulus G" on the strain 
amplitude 70 is also observed at a high oscillation frequency [Fig. [6] (c)]; the weak shear 
thickening occurs at a small strain amplitude 70 < 0.1, and then, the shear thinning occurs at 
a large strain amplitude 70 > 0.1. The storage modulus G' decreases more rapidly than the 
loss modulus G". These features of the shear thinning behaviors can explain the crossover 
behavior of the local storage modulus and loss modulus shown in Fig. [5] (a). At the high 
oscillation frequencies, both local moduli are close to the linear moduli because the local 
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strain amplitude is smaller than about 2 %. In the linear regime, the storage modulus G' is 
larger than the loss modulus G" at high oscillation frequencies. The local storage modulus 
G' deviates from the linear modulus as the oscillation frequency ujq decreases, and the local 
strain 70 exceeds about 2%, while the local loss modulus G" remains close to the linear 
modulus at 70 ~ 2%. Thus, the storage modulus G' crosses over the loss modulus G" at 
a certain oscillation frequency, say a crossover frequency uJq. The loss modulus G" also 
starts to decrease as the oscillation frequency ujq is smaller than the crossover frequency uJq. 
However, because the storage modulus G' decreases more rapidly than the loss modulus G" 
as the strain amplitude 70 increases, the storage modulus G' is smaller than the loss modulus 
G" for oscillation frequencies smaller than the crossover frequency, loq < lOq. 

Figure [7] shows the local Deborah numbers, De" and De^, versus the oscillation frequency 
Wo- The local Deborah numbers, De" and De^, are defined by the products of the oscillation 
frequency ujq and the shear-dependent a relaxation time Ta(7o), De" = ijJoTa{%), and the 
shear-dependent Rouse relaxation time 7-^^(70), De^ = ci;oTr(7o)) respectively. Here, we use 
the fitting functions for the relaxation times and r/j for the simple shear flows obtained 



in Ref. 
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It is seen that the local Deborah numbers De" are equal to unity around the 
crossover frequencies for each local moduli cJq. In the lower figure in Fig. |5] (a), the local 
strain rate 70 increases as the distance from the plate y decreases, while it does not so much 
depend on the oscillation frequency uq but only slightly decreases as the oscillation frequency 
ojq increases. The a relaxation time Tai'j) is the monotonically decreasing function on the 
strain rate 7. Hence, in Fig. [7] (a), the local Deborah number De" decreases as the distance 
from the plate y decreases but does not alter the shape of the curve very much, such that 
the oscillation frequency at which the local Deborah number De" is equal to unity, i.e., the 
crossover frequency cUq, shifts to a higher value as the distance y decreases. As the local 
Deborah number De°' is less than unity, the loss modulus G" also decreases as the oscillation 
frequency Uq decreases as is seen in Fig. |2]for the linear moduli. 

The behaviors of the storage modulus and loss modulus and the amplitude of the local 
strain and strain rate versus oscillation frequency near the plate are shown in Fig. [5] (b). 
The dependency of the local strain and strain rate, 70 and 70, on the oscillation frequency Uq 
near the plate is quite different from those far from the plate. The local strains 70 at y=88 
and 185 slightly increase with the oscillation frequency in the low oscillation frequencies, 
e.g., coq < 0.002, but decrease as the oscillation frequency increases in uq > 0.002, while the 
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local strain at y=10 monotonically increases with the oscillation frequency. The local strain 
rates 70 monotonically increases with the oscillation frequency in the close vicinity of the 
plate, i.e, y=88 and 10, while the strain rate at y=185 shows the non-monotonic dependency 
on the oscillation frequency. The spatial variations of local strains and local strain rates are 
small at low oscillation frequencies and are large at high oscillation frequencies. This feature 
also holds for the local moduli. 

In Fig. [7] (b), we show the local Deborah numbers near the plate. The non-monotonic 
dependency of the local storage modulus G' on the oscillation frequency uq, which can 
be seen in Fig. [S] (b), e.g., G' at y=88 and 185, might be related to the local Deborah 
number De^ in Fig. [7] (b). As can be seen in Fig. E] (b), the local strain rate 70 rapidly 
increases as the oscillation frequency Uq increases at low oscillation frequencies; however, 
at high oscillation frequencies, it shows a different dependency on the oscillation frequency 
Uq. The local strain rate slightly decreases at y=185, does not change much at y=88, 
and monotonically increases at y=10, as the oscillation frequency uq increases. The Rouse 
relaxation time monotonically decreases as the strain rate 70 increases. Hence, the 
local Deborah number De^, De^ = u;oTij(7o), decreases even if the oscillation frequency 
Uq increases at low oscillation frequencies because the Rouse relaxation time r^j rapidly 
decreases. On the contrary, at high oscillation frequencies, De^ increases at y=88 and 185 
as the oscillation frequency Uq increases because tji dose not change much or rather increases 
slightly. The elasticity grows as the local Deborah number De^, but it becomes negligible 
as De^ is less than unity. Thus, the local storage modulus G' at y=88 and 185 in Fig. [5] 
(b) varies according to the local Deborah number De^, but the local storage modulus G' at 
y=10 becomes negligibly small except at low oscillation frequencies. 

In Fig. [HI we also show shear thinning behavior of the local dynamic viscosity r]'{=G" /uq) 
to the local strain rate 70 for various oscillation frequencies, i.e., u)o=1.5 x 10"'^, 6.1 x 10~^, 
1.2 X 10^^, and 2.5 x 10^^. It is seen that, at large strain rates, e.g., 70 > 0.01, the dependence 
of the dynamic viscosity rj' on the oscillation frequency uq is weakened due to which the 
results for different oscillation frequencies coincide with each other for 70 ^ 0.1. A second 
Newtonian regime is also observed at very large strain rates, e.g., 70 > 0.4. We also find that 
the slope of the shear thinning is similar to that observed in steady shear fiows, in which 
the slope is about -0.7—. Thus, the shear thinning behavior similar to that of the steady 
shear flows is observed in the vicinity of the plate at high oscillation frequencies. 



12 



We also show a diagram of the loss tangent tan 5, which is defined as the ratio of the 
viscosity to the elasticity, tan 5 = C'/G", for different oscillation frequencies and local strain 
rates in Fig. [HI In the diagram, the dashed and dotted lines show that the inverses of the 
shear-dependent a and Rouse relaxation times, Taici)^^ and t/j(7)^^, equal the oscillation 
frequency Wq! i-G-; on the dashed and dotted lines, the local Deborah numbers -De" and 
De^ are equal to unity, respectively. The upper side of the diagram indicates the smaller 
distance from the plate while the lower side the larger distance because the local strain rate 
decreases as the distance from the plate increases. No symbols are plotted for large and 
small strain rates at low oscillation frequencies because no data are available given that the 
spatial variations of local strain rates are small at low oscillation frequencies as we have seen 
in Fig. O Below the dashed line, the loss tangent is less than unity, while above the dashed 
line, the loss tangent is larger than unity. Thus, the crossover of the storage modulus and 
loss modulus takes place at the dashed line. Near or above the dotted line, the loss tangent 
is quite large, and the elasticity may be negligible. Hence, the melt behaves as a viscous fluid 
provided that the value of the local strain rate is larger than or close to values that lie on 
the dotted line. Thus, the melt forms three different rheological regimes, i.e., the solid-like 
viscoelastic, liquid-like viscoelastic, and viscous fluid regimes, according to the local strain 
rates and oscillation frequencies. 

D. LAOS analysis 

As the local strain increases near the oscillating plate, the amplitude of higher harmonics 
of the local macroscopic quantities becomes large, and thus, the time evolution of the local 
quantities deforms evidently from that of pure sinusoidal curves. In the present subsection, 
we carry out LAOS analysis^"— to examine the non-linear effects of higher harmonics on 
the local macroscopic quantities. 

Figure [10] shows the power spectra of the local shear stress and strain rate, |5"2.y(a;)p 
and |7(c<;)p, in the rapidly oscillating plates with an oscillation frequency = 0.025. Here, 
Oxyiyj) and 7(0;) represent the Fourier coefficients of the shear stress and strain rates, respec- 
tively. The peaks of the higher harmonics are detected in the odd harmonics, Swq, Scjq, ■ ■ ■ , 
near the oscillating plate but disappear far from the plate. The higher harmonics arise both 
in the local stress and strain rate because the higher harmonics of the local strain rate is 
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induced by the local stress with the higher harmonic contribution and the higher harmonics 
of the local stress is also induced by that local strain rate. 

Figure [TT] shows the spatial variation of the fractional amplitude of the third harmonics 
to the basic oscillation for various oscillation frequencies Uq at Tq = 0.5. The horizontal axis 
is the normalized distance with respect to the thickness of boundary layer shown in Fig. SI 
Even though the local strain monotonically increases while approaching the oscillating plate, 
the fractional amplitude of the third harmonics rather decreases rapidly in the boundary 
layer, i.e., y < lb- This is caused by the boundary effect because the velocity at the oscillating 
plate is purely sinusoidal, and higher harmonics are not allowed in the local strain rate on 
the oscillating plate. Thus, the contribution of the higher harmonics for the local stress 
is depressed in the boundary layer. The fractional amplitude of the third harmonics takes 
the maximum value at the outside of the boundary layer for each oscillation frequency, 
y > lb, and then decreases as the normalized distance increases. We also show the fractional 
amplitude of the third harmonics versus local strain amplitude at different strain amplitudes 
of the system, ro=0.5 and 0.02, for various oscillation frequencies uq in Fig. [121 For a large 
strain amplitude of the system, Fq = 0.5, the fractional amplitude of the third harmonics 
has a maximum value around the position where the local strain amplitude is unity, 70 ~ 1, 
and rapidly decreases in the thickness of boundary layer, y ^ h- On the contrary, for a 
small strain amplitude of the system Fq = 0.02, the maximum occurs not at 70 ^ 1 but 
around the position that the value of the distance from the plate coincides with the value of 
the thickness of boundary layer, y ~ For Fq = 0.02, the position at which the local strain 
amplitude is unity lies at the inside of the boundary layer. Thus, the intrinsic maximum of 
the fractional amplitude of the third harmonics for 70 ~ 1 is suppressed inside the thickness 
of boundary layer, y < lb- 

Figure [13] shows the Lissajous-Bowditch curves of local shear stress a^y versus local strain 
7 at different oscillation frequencies of the plate coq with Fq = 0.5. The dotted line (green) 
indicates a perfect ellipse drawn by the storage and loss modulus calculated from the Fourier 
coefficients of the first harmonics. Thus, the deviation of the Lissajous-Bowditch curve from 
the pure ellipse represents the contribution of the higher harmonics in each figure. At each 
oscillation frequency, the deviation is more evident at some distance away from the plate 
than in close vicinity to the plate. This agrees with the fact that the contribution of the 
higher harmonics is depressed inside the boundary layer and that the fractional amplitude 
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of the third harmonics assumes a maximum at the outside of boundary layer. In Fig. 
[T3l the thicknesses of the boundary layer If, are /fe=223, 45, and 20 for a;o=1.5 x 10^^, 
6.1 X 10^^, and 0.025, respectively. We also show the minimum-strain modulus G^^ and 
large-strain modulus G'^ defined as G\j- = daxy/d'y\^=o and G'^ = cra.j^/7|^=^(,,— respectively. 
The minimum-strain modulus G'j^j yields the elasticity at the point where the change of strain 
rate is zero, d'j/dt = 0, and the large-strain modulus G'^ yields the elasticity at the point 
where the instantaneous strain rate is zero, 7 = 0. Both the minimum-strain G'^ and large- 
strain moduli G^ coincide with the elastic modulus G' in the linear regime, G^j,G^ — )■ G' 
for small 79. It is seen that the large-strain modulus G'^^ is larger than the minimum-strain 
modulus G'nj within a cycle of Lissajous-Bowditch curves at the close vicinity of the plate, 
e.g., y = 10. Thus, intra-cycle stiffening occurs in close vicinity of the plate, although the 
storage modulus G' exhibits shear-thinning behavior for the local strain amplitude 70, i.e., 
inter-cycle softening (see also Fig. |6]). 

Finally, we show the 3D Lissajous-Bowditch curve of the local shear stress a^y, local 
strain 7, and local strain rate 7 near the rapidly oscillating plate. The 2D projections on 
each plane are also shown. The cycle of 7- axy curve is very narrow, and a secondary 
loop^^ is observed at a large instantaneous strain rate. The 7-7 curve also deviates from 
conforming to a perfect ellipse because the higher harmonics are also involved in the local 
strain and strain rate. 

IV. SUMMARY 

We investigated the dynamic rheology of a model polymer melt in non-uniform oscillatory 
shear flows under the transient force between rapidly oscillating plates by using a hybrid 
simulation of the molecular dynamics and computational fluid dynamics. In the quiescent 
state, the melt is in a supercooled state, where the stress relaxation function exhibits a 
stretched exponential form on the time scale of the a relaxation time Tq, and then follows 
the Rouse relaxation function characterized by the Rouse relaxation time tji. [See Fig. [1]] 

In the rapidly oscillating plates, the melt forms a thin boundary layer over the plates 
due to the transient force [See Fig. H] such that the dynamic rheology of the melt spatially 
varies considerably according to the local flow field. At a high oscillation frequency, the 
melt forms three different rheological regimes, i.e., the viscous fluid regime (G' ^ G"), the 
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liquid-like viscoelastic regime (G" < G"), and the solid-like viscoelastic regime {G' > G"), 
over the oscillating plates according to the local Deborah number, while, at a low oscillation 
frequency, the spatial variation of the storage modulus and loss modulus is rather small, 
and the loss modulus is larger than the storage modulus, G' < G", at any distance from the 
plate. 

The dependency of the local moduli on the oscillation frequency, G' and G" vs. uq 
at a fixed position changes according to the distance from the plate. Far from the plate 
[Fig. |5](a)], the local strain 70 decreases as the oscillation frequency ojq increases such that 
the local moduli deviate from the linear moduli larger at a low frequency than at a high 
oscillation frequency. Near the plate [Fig. \5\ (b)], the dependency of the local strain and 
strain rate on the oscillation frequency is quite different from that far from the plate, e.g., 
in close vicinity to the plate, the local strain and strain rate monotonically increase as the 
oscillation frequency. As the local strain rate becomes larger than about 0.01 near the plate, 
the shear thinning of the dynamic viscosity is increased, and the shear thinning behavior 
becomes similar to that observed in steady shear flows [see Fig. [8l] 

The diagram of the loss tangent of the melt for different oscillation frequencies and local 
strain rates is also shown in Fig. [91 It is seen in the diagram that the melt generates different 
rheological regimes according to the oscillation frequency and local strain rates. 

We also investigate the non-linear rheological properties in the LAOS regime in the vicin- 
ity to the oscillating plate. The odd higher harmonics, Suq, Swq, ■ ■ ■ , are detected in the 
power spectra of the local macroscopic quantities near the oscillating plate [see Fig. ITOl ] 
The fractional amplitude of the higher harmonics rapidly decreases inside the boundary layer 
while approaching the oscillating plate, although the local strain increases monotonically. 
[See Fig. [TTl] This is because higher harmonics are not allowed in the velocity on the oscil- 
lating plate due to a non-slip boundary, and thus, the higher harmonics of local quantities 
are suppressed within the boundary layer. 

The Lissajous-Bowditch curve of the local shear stress versus local strain changes the 
shape due to the contribution of the higher harmonics depending on the distance from the 
plate. The elastic modulus exhibits shear-thinning behavior between the different positions 
(inter-thinning) [see Fig. [6], but shear-thickening behavior is observed in a cycle of the 
Lissajous-Bowditch curve at a fixed position in the vicinity of the plate (intra-thickening) 
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[see Fig. [13]. 
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FIG. 1: (Color online) The stress relaxation function G{t) of the model polymer melt in the 
quiescent state (the solid line). The dotted line shows the stretched exponential form Gs{t) in Eq. 
(HI) with Go=5, Ts=90, and c=0.5. The dashed line shows that the Rouse relaxation function GR{t) 
in Eq. ([5]) with r/j=6 x 10^. is the a relaxation time, Ta=3W. 
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FIG. 2: (Color) The storage and loss moduli, G'{lo) (red solid line) and G"{uj) (blue solid line), 
respectively, in the linear response regime. The dashed line shows the slope of 1 as a guide. The 
dotted lines show the shear moduli calculated from the superposition of the stretched exponential 
form Gs{t) and the Rouse relaxation function Gji{t), Gs{t) + Gii{t). The diamonds aligned on 
the upper horizontal axis indicate the values of the oscillation frequencies to be used in the hybrid 
simulations. 
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FIG. 3: Schematics of the geometry of the problem, the mesh system, and the time-evolution 
scheme. 
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(a) velocity profile (b) 4 vs. ojq 

FIG. 4: (Color online) The snapshots of the velocity profile for a;o=0.0015 and ro=0.5 (a) and the 
thickness of boundary layer If, for various oscillation frequencies ojq (b). In figure (b), the solid 
line, dashed line with squares, and dashed line with triangles show the results for the linear regime 
(Fq ^ 1), Fo=0.5, and Fo=0.02, respectively. The dotted line in figure (b) shows the slope of 
for the Newtonian fluid. 
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(a) Far from the plate (b) Near the plate 



FIG. 5: (Color) The local dynamic rheology, G', G" vs. wq, and the amplitude of local strain and 
that of the local strain rate, 70 and 70, respectively, versus the oscillation frequency ujq far from the 
plate (a), e.g., y=1309 (Blue), y=771 (Green) and 381 (Red), and near the plate (b), e.g., y = 185 
(Blue), y=88 (Green) and 10 (Red). The black lines show the linear moduli depicted in Fig. [21 
The open circles in the upper figure of (a) show the crossover points of G' and G" . The dotted 
lines in the lower figure of (a) show 70=0. 02. 
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(a) Wo = 6.1 X 10 



(b) uo = 1.2 X 10" 
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(c) Luo = 2.5 X 10-2 



FIG. 6: The local storage and loss moduli, G' and G" , vs. the local strain amplitude 70 for 
uo=6.1 X 10-3 (a), 1.2 X Q-^ (b), and 2.5 x lO'^ (c). 
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FIG. 7: (Color) The local Deborah numbers defined as De°'=ujo'Ta{'jo) and De^=ujoTji{'y(j) far from 
the plate (a) and near the plate (b). Here, Ta^R^jo) are the a and Rouse relaxation times depending 
on the local strain rate 79. We use the formulas for Ta^ni'y) that are obtained by the molecular 
dynamic simulations of the model polymer melt in steady shear flows2S. The solid lines show De" 
and the dashed lines De^. The black lines show the Deborah numbers defined via the oscillation 
frequency ojq and the relaxation times of the melt in the quiescent state Ta^miO)- The downward 
arrows show the crossover frequencies cjq at each distance from the plate, which are also shown by 
open circles in Fig. [5] (a). 
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FIG. 8: The dynamic viscosity r/' vs. the ampUtude of the strain rate 70 for ljJq=3.1 x 10 
1.2 X 10-2, and 2.5 x 10"^. Here, rj' = G"/oj and 7=tJo7o- 
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FIG. 9: (Color) The amplitude of the loss tangent tan 5, which is defined as tan 5 = G"/G', for 
different oscillation frequencies and local strain rates. Here, G' and G" are averages of the local 
storage modulus and loss modulus provided that the local strain rate is in each interval with 0.5 
logarithmic scale centered at the value of the local strain rate of each plot, respectively. On the 
dashed and dotted lines, the inverses of the shear-dependent a and Rouse relaxation times, Tq,(7) 
and T/j(7), equal the oscillation frequency, respectively. 
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FIG. 10: The power spectra of the local shear stress and strain rate, Idxyiio)]"^ (the upper figures) 
and |7(cj)|^ (the lower figures), at difi^erent positions for ujq = 0.025 and Fq = 0.5. 
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FIG. 11: (Color online) The spatial variation of the fractional amplitude of the third harmonics 
to the first harmonics of the local shear stress, \^xyi'i^o)\'^ /\(^xyi^o)\'^ , for various oscillation fre- 
quencies of the plate coq at Fq = 0.5. The horizontal axis shows the distance normalized to the 
thickness of the boundary layer at each oscillation frequency ojq, shown in Fig. HI 
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FIG. 12: (Color online) The fractional amplitude of the third harmonics to the first harmonics 
of the local shear stress, \d'xy{3u)o)\'^ /\^xy{^o)\'^ , versus the amplitude of local strain, 70, at the 
different amplitudes of strain of the system ro=0.5 and 0.02 for various oscillation frequencies uq. 
The shape of the symbol represents the oscillation frequency ojq, which is the same as in Fig. [TTl 
The downward and upward arrows show the position oi y = Ij, for each oscillation frequency of the 
plate at Fq = 0.5 and Fq = 0.02, respectively. 
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FIG. 13: (Color) The Lissajous-Bowditch curves of the local shear stress CTxy vs. the local strain 7 
at different oscillation frequencies of the plate loq with Fq = 0.5. The dashed line (blue) shows the 
minimum-strain modulus (or tangent modulus) G'j^, G'j^j = daxy/d'y\-y=o, and the dash-dotted line 
(red) shows the large-strain modulus (or scant modulus) G'^, G'^ = axy/j\'y=-yg- The dotted curve 



(green) shows the pure ellipse formed by the storage and loss modulus for the first harmonics. 



31 



ro=0.5, wo=3.1xl0-^\ y=29 



■"■^io' ' ' '-I'o' ' ' ' 6 ' ' ' Vo' ' ' io 

7 



FIG. 14: (Color) 3D Lissajous-Bowditch curves of the local shear stress axy, local strain 7, and 
local strain rate 7 near the rapidly oscillating plate. 2D projections onto the planes of the shear 
stress axy and strain 7, stress axy and strain rate 7, and strain rate 7 and strain 7 are also shown 
on the right side. 
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